library(tidyverse)
library(tidytext)
library(lubridate)

source("prep.R")

df <- read_csv("outputs/facts.csv") %>%
  mutate(year = year(date))

ingredients <- read_csv("outputs/ingredients.csv")

Yearly numbers

Let’s see how the number of recipes have been by year.

df %>%
  count(year) %>%
  arrange(year)

df %>%
  count(year) %>%
  ggplot(aes(year, n)) +
  geom_line()

df %>%
  count(year) %>%
  mutate(n = cumsum(n)) %>%
  ggplot(aes(year, n)) +
  geom_line()

There weren’t too many before 2004, a huge spike in 2004, and then a slow taper afterwards. They probably did some backend retooling which re-dated all the articles before 2004 to 2004. It is strange that the number of articles has been steadily decreasing.

Let’s also look at the ratings. Based on what we just saw, we will focus only after 2004.

df %>%
  filter(year > 2004) %>%
  group_by(year) %>%
  summarize(rating = mean(rating, na.rm = TRUE)) %>%
  ggplot(aes(year, rating)) +
  geom_line()

Well the average rating has been going up, so it seems like Epicurious has been focusing on quality over quantity.

There was a serious spike in 2021, which is interesting. The data was scraped in November, so there was enough data that it probably wasn’t just chance. Perhaps this was pandemic related?

Authors

Who are the most prolific authors?

unnest_authors <- df %>%
  unnest_tokens(author, author, token = str_split, pattern = ";")

unnest_authors %>%
  filter(!is.na(author)) %>%
  count(author) %>%
  arrange(-n)

Bunch of Bon Appetit people, not surprising given that Conde Nast owns both. Let’s look at when they were prolific.

df %>%
  filter(!is.na(author)) %>%
  add_count(author) %>%
  filter(n >= 169) %>% # get top 10
  count(author, year) %>%
  ggplot(aes(year, n)) +
  geom_line() +
  facet_wrap(~author, scales = "free_y")

This all makes sense. BA published as a team until around 2016, which is also when some of the BA people started showing up. They were probably pushing the personalities more. Also note that a lot of the BA people stopped showing up around last year, which was when the blow up happened.

We also look at the most prolific authors by year.

unnest_authors %>%
  filter(year > 2004, !is.na(author)) %>%
  count(year, author) %>%
  group_by(year) %>%
  slice_max(n, n = 5)

Ingredients

We now analyze the ingredients of recipes. We first expand the table by each ingredient in the recipe.

ingr <- df %>%
  left_join(ingredients, by = "id") %>%
  prep(to_mat = FALSE)

ingr

Base is the cleaned ingredient name. Let’s count how times each ingredients have shown up. What are the most popular ingredients?

ingr %>%
  count(base) %>%
  arrange(-n) %>%
  head(100)

What is the distribution of the number of ingredients per recipe like?

ingr %>%
  count(id) %>%
  ggplot(aes(n)) +
  geom_histogram()

`stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

We will look at how ingredients change under two divisions: time and rating. Let’s first look at time. We will split time into two parts: early <= 2004 (due to the timestamp in the data) and late >= 2014 (when Epicurious and BA were combined in a digital platform). We want to calculate log odds, so we will ignore the inbetween.

log_odds <- function(df, field, a, b) {
  df %>%
    count(base, {{ field }}) %>%
    group_by({{ field }}) %>%
    mutate(prob = (n + 1) / (sum(n) + 1)) %>%
    pivot_wider(base, names_from = {{ field }}, values_from = prob, values_fill = 0) %>%
    filter({{ a }} > 0, {{ b }} > 0) %>%
    mutate(log_odds = log({{ a }}, base = 2) - log({{ b }}, base = 2))
}

ingr %>%
  mutate(period = case_when(
    year <= 2004 ~ "Early",
    year >= 2014 ~ "Late",
    TRUE ~ "Middle"
  )) %>%
  filter(period != "Middle") %>%
  log_odds(period, Late, Early) %>%
  {
    \(x) bind_rows(slice_max(., log_odds, n = 10), slice_min(., log_odds, n = 10))
  }() %>%
  mutate(base = fct_reorder(base, log_odds)) %>%
  ggplot(aes(log_odds, base)) +
  geom_col(aes(fill = log_odds > 0)) +
  theme(legend.position = "none") +
  labs(x = "Log odds of appearing in later recipes", y = "Ingredient")

We also look at the difference between the top and bottom quartiles. Recipes of the bottom quartile tend to feature basic ingredients like pastas and legumes whereas the top rated recipes have more exotic ingredients like ramps, cotija, and grappa.

quantile(df$rating, probs = c(0.25, 0.75), na.rm = TRUE)

 25%  75% 
3.00 3.66 

ingr %>%
  mutate(group = case_when(
    rating <= 3 ~ "Lo",
    rating >= 3.66 ~ "Hi",
    TRUE ~ "Middle"
  )) %>%
  filter(group != "Middle") %>%
  log_odds(group, Hi, Lo) %>%
  {
    \(x) bind_rows(slice_max(., log_odds, n = 10), slice_min(., log_odds, n = 10))
  }() %>%
  mutate(base = fct_reorder(base, log_odds)) %>%
  ggplot(aes(log_odds, base)) +
  geom_col(aes(fill = log_odds > 0)) +
  theme(legend.position = "none") +
  labs(x = "Log odds of appearing in highly rated recipes", y = "Ingredient")

Similarities

We now look at the similarity matrix that we will use in the recommender.

source("ingredients-rec/ingr.rec.R")

sim <- ingr %>%
  prep(min_num = 10) %>%
  create_sim()

We will pivot the similarity table into long form.

sims <- sim$sim %>%
  as_tibble() %>%
  mutate(ingr1 = colnames(.)) %>%
  pivot_longer(!ingr1, names_to = "ingr2", values_to = "val") %>%
  filter(ingr1 < ingr2)

sims

Let’s look at a histogram of the similarities.

sims %>%
  ggplot(aes(val)) +
  stat_ecdf()

Over 80% of pairs have no similarity value and most values are under 0.2.

We now look at a graph where we connect ingredients with particularly high similarty.

library(igraph)
library(ggraph)

set.seed(2022) # reproducibility

sims %>%
  filter(val > 0.22) %>%
  graph_from_data_frame() %>%
  ggraph(layout = "fr") +
  geom_edge_link() +
  geom_node_point() +
  geom_node_text(aes(label = name), size = 3, vjust = 1, hjust = 1)

We see some natural connected components like baking, asian stir fry, cocktail, and gluten free baking.

Clustering

Finally, we will do clustering of the item vectors. We first create the item vector matrix where the vectors are all unit. This is code from the create_sim function.

mat <- df %>%
  left_join(ingredients, by = "id") %>%
  prep(min_num = 5)

mat <- t(mat / sqrt(rowSums(mat)))
mat <- mat / sqrt(rowSums(mat * mat))

We now run kmeans to cluster the vectors. Note that Euclidean distance between unit vectors corresponds well with the cosine of the angle by law of cosines.

set.seed(2021) # reproducibility

res <- kmeans(mat, centers = 10, nstart = 10)

tibble(
  ingredient = rownames(mat),
  cluster = res$cluster
) %>%
  group_by(cluster) %>%
  slice_head(n = 5) %>%
  print(n = 50)

# A tibble: 50 × 2
# Groups:   cluster [10]
   ingredient              cluster
   <chr>                     <int>
 1 table salt                    1
 2 fine sea salt                 1
 3 avocado oil                   1
 4 medjool date                  1
 5 yeast                         1
 6 chamomile tea bags            2
 7 fruit                         2
 8 creme de cassis               2
 9 grapes                        2
10 mint leaves                   2
11 button mushrooms              3
12 oyster mushrooms              3
13 chiles de arbol               3
14 sherry wine vinegar           3
15 horseradish                   3
16 black peppercorn              4
17 lemongrass                    4
18 plain whole milk yogurt       4
19 xanthan gum                   4
20 coriander                     4
21 angostura bitters             5
22 orange twist                  5
23 light rum                     5
24 vodka                         5
25 ice cube                      5
26 white onion                   6
27 sour cream                    6
28 cilantro leaves               6
29 lime juice                    6
30 avocado                       6
31 honey                         7
32 molasses                      7
33 maple syrup                   7
34 almond                        7
35 chocolate                     7
36 red bell pepper               8
37 olive oil                     8
38 garlic                        8
39 bay leave                     8
40 worcestershire sauce          8
41 salt                          9
42 unsalted butter               9
43 egg                           9
44 sugar                         9
45 vanilla extract               9
46 shiitake mushroom            10
47 vegetable oil                10
48 soy sauce                    10
49 asian sesame oil             10
50 hoisin sauce                 10

Some reasonable clusters. Cluster 5 is obviously cocktails, 7 is probably sweet baked goods, 9 is more generic baked goods, and 10 is Asian cuisine. The other clusters also have themes, but less obvious to interpretation. It seems there is definitely information in the distance between item vectors.